Signatures of Magnetized Large Scale Structure in Ultra-High Energy Cosmic Rays 
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We investigate the impact of a structured universe in the multi-pole moments, auto-correlation 
function, and cluster statistics of cosmic rays above 10 19 eV. We compare structured and uniform 
source distributions with and without magnetic fields obtained from a cosmological simulation. We 
find that current data marginally favor structured source distributions and magnetic fields reaching 
a few /iG in galaxy clusters but below 0.1 fiG in our local extragalactic neighborhood. A pronounced 
GZK cutoff is also predicted in this scenario. Future experiments will make the degree scale auto- 
correlation function a sensitive probe of micro Gauss fields surrounding the sources. 
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PACS numbers: 98.70.Sa, 13.85.Tp, 98.65.Dx, 98.54.Cm 

Introduction: The origin and nature of ultra high en- 
ergy cosmic rays (UHECRs) is one of the most puzzling 
problems of high energy astrophysics . Beside the ac- 
tual mechanism responsible for the acceleration of the 
particles to energies as high as a few x 10 20 eV and more, 
propagation of UHECRs in intergalactic space is a deli- 
cate and critical part of the full investigation. It was soon 
recognized that energy losses dramatically limit the pos- 
sible distance traveled by these particles; in particular, 
interactions with cosmic microwave background photons 
lead to catastrophic photo-meson production at energies 
above ~5x 10 19 eV, the GZK effect 0. Charged nucle- 
ons are also susceptible to deflections induced by extra- 
galactic magnetic fields (EGMF). The importance and 
the subtlety of the effects produced by this process have 
been emphasized by a number of authors 0, 0, 0, @ 
However, up to date, the properties of EGMF are only 
very poorly constrained; therefore, their impact on the 
propagation of UHECRs is still quite uncertain and a 
number of different scenarios, compatible with the avail- 
able observational constraints, must be investigated. 

Along this line of reasoning, recently we have investi- 
gated UHECR propagation in a highly structured EGMF 
obtained from a simulation of large scale structure for- 
mation 0. The EGMF was evolved according to the 
induction equation. Its seeds were generated at shocks 
in accord to the Biermann battery mechanism 8]. Thus, 
it was amplified in different parts of the universe accord- 
ing to the velocity field provided by the simulated gas 
component. Lorentz forces were neglected as EGMF are 
known to be dynamically unimportant in most of the uni- 
verse volume. The resulting field strength was normal- 
ized so as to reproduce observed rotation measure values 
in galaxy clusters The EGMF are significant only 
within filaments and groups/clusters of galaxies. The 
relative strength there and the resulting overall topol- 
ogy depend solely on hydrodynamics which is basically 
driven by gravitational instability. About 90 percent of 
the volume is filled with negligible fields. This is a reflec- 
tion of the adopted Biermann field generation scheme, as 
opposed to the case in which the initial magnetic field is 



set uniform over the whole volume. 

Here we extend our previous study Q to include 
sources at cosmological distances, model more realisti- 
cally the UHECR source distributions and compare ex- 
plicitly with the case of no EGMF. Our favored scenario 
is in considerably improved agreement with UHECR data 
and reveals signatures of the EGMF to which current cos- 
mic ray data start to become sensitive. 

Numerical Techniques: The simulation of large scale 
structure formation was carried out within a volume of 
size 50ft. _1 Mpc, on a side, where h = Ho/(100 km s _1 
Mpc -1 ) = 0.67, using a comoving grid of 512 3 zones and 
256 3 dark matter particles ^J. Two observer positions 
are compared as in Ref. Q: The first in a filament-like 
structure with EGMF ~ 0.1 /iG and the second at the 
border of a small void with EGMF ~ 10" 11 G. A rela- 
tively large structure about 17 Mpc away from the weak 
field observer is identified for calculation purposes as the 
Virgo cluster. We orient our terrestrial coordinate sys- 
tem so that this cluster is close to the equatorial plane. 

For a given number density of UHECR sources, n s , 
we explore both the case in which their spatial distri- 
bution is either proportional to the local baryon den- 
sity, as in Ref. [jj, or completely uniform. Further- 
more, we assume that these sources are distributed (1) 
with respect to their (dimensionless) UHECR emission 
power, Qi, as dn s /dQi oc Q~ 2 ' 2 with 1 < Qi < 100, 
as motivated by the luminosity function of the EGRET 
7— ray blazars in the power range 10 46 ergs _1 < Qi < 
10 48 ergs _1 And (2) with respect to the spectral in- 
dex cti of the emitted power-law distributions of UCHRs 
as dn s /doti =const for a — 0.1 < on < a + 0.1, where 
a (~ 2.5) is a free parameter to be determined through a 
best fit analysis. Each source accelerates up to 10 21 eV, 
and its average UHECR emission characteristics do not 
change significantly on the time scale of UHECR propa- 
gation. Injected power and spectral index are then fit to 
reproduce the observed spectrum. 

Taking advantage of the periodic boundary conditions, 
periodic images of the simulation box are added until the 
linear size of the enclosed volume is larger than the en- 
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ergy loss length of nucleons above 10 19 eV, ~ 1 Gpc. As 
in previous work 0, , particles are propagated taking 
into account Lorentz forces due to the EGMF and en- 
ergy losses. In this respect, pion production is treated 
stochastically, whereas pair production as a continuous 
energy loss process. Trajectories connecting sources and 
observers in different copies of the simulation box are 
taken into account. 

An event was registered and its arrival directions and 
energies recorded each time the trajectory of the propa- 
gating particle crossed a sphere of radius 1.5 Mpc around 
one of the observers. For each configuration 10 4 events 
(a realization) where collected to construct simulated sky 
distributions. These are multiplied with the solid-angle 
dependent exposure function and folded over the angular 
resolution characteristic of a given experiment. Finally 
from each simulated sky distribution typically 100 mock 
data sets, consisting of iV bs observed events, were ran- 
domly extracted. The angular resolution is ~ 1.6° above 
4 x 10 19 eV and ~ 2.5° above 10 19 eV for the AGASA 
experiment []3 and ~ 10° for the SUGAR data 0. 
For the exposure function we use the parameterization 
of Ref . with parameters as in Ref . [l5| . Energy reso- 
lution effects, supposed to be of order AE/E ~ 30%, are 
also taken into account. 

For each data set (both simulated and real) we 
obtained estimates for the spherical harmonic coef- 
ficients C(Z), the autocorrelation function N(9), and 
the number of multiplets M(n) of n events within an 
angle 9 m via the following estimator functions C(Z) = 

(21 - 



where u>i is the total experimental exposure at 
arrival direction u l , and Yi m (u l ) is the real- 
valued spherical harmonics function taken at 
direction u l 7]; N(0) is defined as N(0) oc 

{1 ^ % is m same bin as 9, otherwise} /S{9), 
where S(9) is the solid angle size of the corresponding 
bin, and normalized to one for an isotropic distribution. 

The mock data sets in the various realizations yield the 
statistical distributions of C(l), N(9), and M(n). We de- 
fine the average over all mock data sets and realizations 
and two errors. A smaller, statistical error, i.e. the fluc- 
tuations due to the finite number N a b s of observed events, 
averaged over all realizations. And a larger, "total error" , 
i.e. the statistical error plus the cosmic variance due to 
the variation in source characteristics. Typically, for each 
scenario we consider 10 realizations for the source posi- 
tions for each of which we chose 50 different realizations 
for the power Qi, and injection spectral index o^, sam- 
pled from the distributions discussed above. 

Given a set of observed and simulated events, we define 

X4 = Z). [(Si,data - Si.simu) /A<Sj )S i m u] , where Si stands 
for Ci, N(9), or M(n), Sj.data refers to 5» obtained from 
the real data, and Si iS - lmu and AS , i iS ; mu are the average 
and standard deviations of the simulated data sets. This 
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measure of deviation from the average prediction can be 
used to obtain an overall likelihood for the consistency 
of a given theoretical model with an observed data set 
by counting the fraction of simulated data sets with \4 
larger than the one for the real data. For further details 
see Refs. HUES [3. 

Constraints from Existing Data and Predictions for 
Future Experiments: The scenarios studied and results 
are summarized in Tab. [I] UHECR sources whose num- 
ber density is given in column II are distributed cither 
proportionally to the simulated baryonic density or uni- 
formly ( "yes" or "no" respectively in column III) . Finally, 
the EGMF is cither taken from the simulation with local 
value as indicated in column IV, or completely neglected 
("no EGMF"). The number of realizations of source po- 
sitions was 10 except for scenario 1 for which it is 7. 
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FIG. 1: Angular power spectrum C(l) as a function of multi- 
pole I, obtained for the AGASA+SUGAR exposure function, 
for iV bs = 99 events observed above 40 EeV. We show the 
realization averages (diamonds), statistical (left) and total 
(right) error bars, respectively, predicted by the model. The 
red (thick) diamonds and red (thick, outer) error bars repre- 
sent scenario 6, whereas the black (thin) diamonds and black 
(thin, inner) error bars represent scenario 1. The histogram 
corresponds to the AGASA+SUGAR data. 

It was noted in Ref. that the AGASA and SUGAR 
experiments had comparable exposure in the northern 
and southern hemisphere, respectively. Above 4 x 10 19 eV 
they combined 50 events observed by AGASA (exclud- 
ing 7 events observed by Akeno) and 49 events seen 
by SUGAR. While SUGAR'S angular resolution is much 
worse than for AGASA and in general prevents a com- 
bination of the two data sets, multi-poles I < 10 are not 
sensitive to scales < 10°. To take advantage of more 
statistics and more complete sky coverage, we therefore 
use the combined AGASA+SUGAR data set of 99 events 
above 4 x 10 19 eV when computing multi-poles I < 10. 
The results for scenarios 1 and 6 in Tab. [fl are compared 
in Fig.n In contrast, since the auto-correlation function 
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TABLE I: List of simulated scenarios. The columns contain the number assigned to the scenario, the source density, whether 
the sources are distributed as the baryon density in the simulation box or uniformly (yes/no), the magnetic field strength at 
the observer location, the best fit power law index in the injection spectrum E~ a , and the overall likelihoods of fits to the 
data. The first four likelihoods are for the multi-poles above the energy indicated as superscript in EeV and over the range 
of I indicated as subscript. Above 40 EeV the comparison was made with iVobs = 99 AGASA+SUGAR events, whereas above 
10 EeV comparison with an isotropic distribution of 1500 events was made, see text. The last two likelihoods are for the 
auto-correlation for 6 < 20°, and the clustering within 2.5° up to multiplicity 10, respectively. 
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and clustering are sensitive to small scales, when com- 
paring with data above 4 x 10 19 eV we will only use the 
reported 57 AGASA+Akeno events. 

No signs of anisotropy were found by either the 
AGASA or the SUGAR experiment down to 10 19 eV. We 
therefore also compute predictions for the multi-poles for 
~ 1500 events observed above 10 19 eV, using the com- 
bined exposure of AGASA+SUGAR. This is compared 
with an isotropic distribution in Fig. |3 for scenarios 1 
and 6. The likelihoods for consistency with isotropy are 
also summarized in Tab.[I] 

N oba ( 10. < E/EeV < 1000.)= 1500 



0.0010 



0.0001 



1 1* 



atively strong fields surrounding the observer suppress 
the flux from the far sources. The observed flux is there- 
fore more strongly dominated by the nearest few sources 
and thus more anisotropic. For observers immersed in 
weak fields, a relatively strong GZK cutoff is expected, 
inconsistent with AGASA data, see Fig. |3 and in con- 
trast to the disfavored strong field observer. 
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FIG. 2: As Fig. but for the angular power spectrum ob- 
tained for the combined AGASA+SUGAR exposure function, 
for iVobs = 1500 events observed above 10 EeV. The straight 
line is the analytical prediction, Ci ~ (4-7T./V) -1 , for the aver- 
age multi-poles for complete isotropy. 

Comparison in Tab. [I] of multi-poles, auto-correlation 
function and clustering above 4 x 10 19 eV, as well as ob- 
served isotropy at 10 19 eV, see Fig. favor an observer 
immersed in relatively weak fields, especially scenario 6. 
For the multi-poles describing the large scale distribu- 
tion this can qualitatively be understood as follows: Rel- 



FIG. 3: Predicted spectrum observable by AGASA for sce- 
nario 6, for which multi-poles were shown in Fig. as com- 
pared to the AGASA (dots) and HiRes-I [T^j (stars) data. The 
histogram marks the average and the two lines above and be- 
low the 1-sigma fluctuations over the simulated realizations. 
The solid straight line marks the injection spectrum. 

Scenarios 3 and 4 are somewhat disfavored by auto- 
correlation and clustering: Structured sources produce 
more clustering in the absence of magnetic fields. While 
this effect is only marginal in the current data set, it 
will provide a solid criterion for discriminating between 
the two scenarios in future experiments: Fig. 01 shows 
that the auto-correlation function is sensitive to the mag- 
netic fields around the sources. Scenarios with no signifi- 
cant magnetic fields predict a stronger auto-correlation at 
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small angles, independent of whether or not the sources 
are structured. This is because the images of sources im- 
mersed in considerable magnetic fields are smeared out, 
which also smears out the auto-correlation function over 
a few degrees. As shown in Fig. 0] this effect will be- 
come quite significant once ~ 1500 events are observed 
above 4 x 10 19 eV, for example, by the Pierre Auger 
project 0|- ^ a ^ so appears robust against variation of 
other parameters such as source density and distribu- 
tion: For example, whereas N(l°) ~ 10.1 ± 1.6 for sce- 
nario 6, all other scenarios without magnetic field pre- 
dict ./V(l°) > 16 ±3. Finally, all scenarios considered 
predict significant anisotropy above 4 x 10 19 eV for expo- 
sures reached by the Pierre Auger project. 

N obs ( 40. < E/EeV < 1000,)= 1500 
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FIG. 4: The auto-correlation function N(6) as a function of 
angular distance 9 for iV b s = 1500 events observed above 
40 EeV in 1° bins. The key for the model predictions is as 
in Figs. Q and |3J now comparing scenario 6 (lower set) and 4 
(upper set) instead, which only differ in the EGMF. N(6) = 1 
corresponds to isotropy. 

Conclusions: Our results in Ref. showed that, in the 
case of structured extragalactic magnetic fields, the data 
were marginally consistent with an observer immersed 
in ~ 0.1 /iG field, although the isotropy of UHECRs at 
10 19 eV could not be explained by a local source dis- 
tribution alone. In the present paper we took into ac- 
count sources up to cosmological distances and adopted 
more realistic distributions for the emission properties of 
UHECR sources. This considerably improved the quality 
of our fits at all energies. We conclude that the combina- 
tion of statistical quantities characterizing the compar- 
ison between observed and simulated data moderately 
favors a scenario in which (i) UHECR sources have a 
density n s ~ 10~ 5 Mpc~ 3 and follow the matter distri- 
bution (ii) magnetic fields are relatively pervasive within 
the large scale structure, including filaments, and with a 
strength of order of a /iG in galaxy clusters (iii) the lo- 
cal extragalactic environment is characterized by a weak 
magnetic field below 0.1 /xG. Thus, current UHECR data 



already allow us to probe the large scale distribution 
of matter and magnetic fields and seem to marginally 
confirm the quite realistic paradigm summarized above! 
This will increasingly be true for future detectors such as 
the Pierre Auger [lg and EUSO projects. In par- 
ticular, the degree-scale auto-correlation functions above 
~ 4 x 10 19 eV can serve as a discriminator between mag- 
netized and unmagnetized sources: N(l°) < 10 points 
to strong magnetization, rather independently of other 
source parameters. The best fit scenario also predicts a 
pronounced GZK cutoff. 
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